open(FH,"/Volumes/HD_2/Documents/workspace/data/BsiSNP_comparison/Normal_SNP.filterCNV.txt.minReadCov_1.snpcaller.test.txt");
my @in=<FH>;
chomp(@in);
close(FH);
my $count = 0;
my $count_het = 0;
my $count_array_not_C = 0;
my $count_hom = 0;
my $fp = 0;
my $fn = 0;
my $total_count=0;
my $ag_count=0;
my $cg_count = 0;
my $count_array_hom_C = 0;
#open(FH,"/Users/yaping/Documents/workspace/gatk_read/Sting/snps.raw.variantonly.vcf");
open(FH,"/Volumes/HD_2/Documents/workspace/snps.raw.performance.newGATK.test.vcf");

my @vcf=<FH>;
chomp(@vcf);
close(FH);


my %vcf_hash = ();
foreach my $line(@vcf){
	next if $line=~/^\#/;
	my @splitvcf=split "\t",$line;
	my $key = $splitvcf[0]."~".$splitvcf[1];
	#next if $splitvcf[5] <=10;
	$vcf_hash{$key} = $line;
}
foreach my $line(@in){
	my @splitin=split "\t",$line;
	if(($splitin[3] eq $splitin[4]) and ($splitin[5] eq "C" or $splitin[5] eq "G")){
				$count_array_hom_C++;
			}
			else{
				$count_array_not_C++;
			}
	my $key = "chr".$splitin[1]."~".$splitin[2];
	if( exists $vcf_hash{$key} ){
			#	$count++;
				my @splitvcf=split "\t",$vcf_hash{$key};
				if(($splitvcf[6] eq "PASS") and (($splitvcf[9]=~/0\/0/ and ($splitvcf[3] eq "C" or $splitvcf[3] eq "G" )) or ($splitvcf[9]=~/1\/1/ and ($splitvcf[4] eq "C" or $splitvcf[4] =~ /\bC\,/ or $splitvcf[4] eq "G" or $splitvcf[4] =~ /\bG\,/)) or ($splitvcf[9]=~/2\/2/ and ($splitvcf[4] =~ /\,C/ or $splitvcf[4] =~ /\,G/)))){
					if(($splitin[3] ne $splitin[4]) or ($splitin[5] ne "C" and $splitin[5] ne "G" and $splitin[5] ne "I")){
						$fp++;
						print("fp:\t$splitin[0]\t$splitin[1]\t$splitin[2]\n");
						#	if(($splitin[5] ne "C")){
							#	print("fp:\t$splitin[0]\t$splitin[1]\t$splitin[2]\n");
						#	}
					
					}
					else{
						$true_pos++;
					}
				}
				else{

						if(($splitin[3] eq $splitin[4]) and ($splitin[5] eq "C" or $splitin[5] eq "G")){
							$fn++;
							
					#		print("false neg:\t$splitin[0]\t$splitin[1]\t$splitin[2]\n");
							
						}
						else{
							$true_neg++;
						}


				}
			}
			else{
				if(($splitin[3] eq $splitin[4]) and ($splitin[5] eq "C" or $splitin[5] eq "G")){
					$fn++;
				
					
					print("non-found-false neg:\t$splitin[0]\t$splitin[1]\t$splitin[2]\n");
				}
				else{
					$true_neg++;
				}
				#else{
				#	$fn++;
				#}
				#$count_hom++;
				##print("$splitin[2]\n");
		}
	$total_count++;
}
my $fp_rate = $fp/$count_array_not_C;
my $fn_rate = $fn/$count_array_hom_C;
print("count is: $count\n");
print("total count is: $total_count\n");
print("het count is: $count_het\n");
print("array cytosine count is: $count_array_hom_C\n");
print("array not cytosine count is: $count_array_not_C\n");
print("false pos num is: $fp\t fp rate is: $fp_rate\n");
print("false neg num is: $fn\t fn rate is: $fn_rate\n");
print("ag count in fn is: $ag_count\n");
print("cpg count is: $cg_count\n");
